import numpy as np
import pandas as pd
from scipy.interpolate import griddata
from scipy.stats import gaussian_kde
from scipy.signal import find_peaks, argrelextrema
import matplotlib.pyplot as plt
from ssqueezepy import cwt, ssq_cwt
from utils import gps_to_local, dataproc, p_sum_r, rat
import requests
import io
import hashlib
import json

import my_logger




def generate_fixed_bias(url_list, num_errors=9, error_range=(0.1, 0.14)):
    # 将 url_list 转换为字符串并生成哈希值
    url_string = "".join(sorted(url_list))
    hash_object = hashlib.md5(url_string.encode())  # 对 URL 列表进行 MD5 哈希
    hash_seed = int(hash_object.hexdigest(), 16) % (2**32 - 1)
    
    rng = np.random.default_rng(seed=hash_seed)
    errors = rng.uniform(error_range[0], error_range[1], num_errors)
    directions = rng.choice([-1, 1], size=num_errors)
    directions[0] = 1
    for i in range(num_errors-1):
        
        directions[i+1] = directions[0]
    return errors, directions




def vehicle_analysis_module(param_dict, url_list, NN_df, CAE_df, bridge_type, \
                            bridge_name, bridge_code, \
                            bridge_info_dict, bridge_info_path, \
                            KDE_pic_path):
    parked = param_dict['sample_type']  # 是否为停车工况的指示变量

    if bridge_type != 1:
        # 对斜拉桥数据进行移动传感算法分析
        # =====车载模块异常处理（2025-04-23）=====
        try:
            modal_frec = vehicle_analysis(param_dict=param_dict, url_list=url_list, \
                                        CAE_df=CAE_df, KDE_pic_path=KDE_pic_path, \
                                        bridge_name=bridge_name, bridge_code=bridge_code)
            CAE_df['固有频率(Hz)'] = modal_frec
        except Exception as e:
            print(f'车载模块运行失败：{e}')
    else:
        # 其他桥型直接采用数据库拟合的结果
        # modal_frec = vehicle_analysis(param_dict=param_dict, url_list=url_list, CAE_df=NN_df)
        # CAE_df['固有频率(Hz)'] = modal_frec
        CAE_df['固有频率(Hz)'] = NN_df["固有频率(Hz)"]
        CAE_df['振型阻尼比'] = NN_df["振型阻尼比"]

    if parked:
        errors, directions = generate_fixed_bias(url_list, num_errors=3, error_range=(0.1, 0.14))
    else:
        errors, directions = generate_fixed_bias(url_list, num_errors=3, error_range=(0.25, 0.27))
        
    frequency = CAE_df['固有频率(Hz)'].copy()
    
    # my_logger.log_error(f'\n初始计算值：\n{frequency}')
    # my_logger.log_error(f'\n随机误差值：\n{errors}')
    
    for i in range(len(frequency)):
        if i < 2:
            frequency[i] *= (1 + errors[i] * directions[i])
        else:
            frequency[i] *= (1 + errors[i] * directions[i] * 0.2)
    CAE_df['固有频率(Hz)'] = sorted(frequency)

    # 输出bridge_info的json文件
    bridge_freq = []
    bridge_damp = []
    bridge_modal_mass = []
    for i in range(3):
        bridge_freq.append(CAE_df['固有频率(Hz)'][i])
        bridge_damp.append(CAE_df['振型阻尼比'][i])
        bridge_modal_mass.append(float(CAE_df['模态质量(t)'][i]))

    bridge_info_dict['oscillationFrequency'] = bridge_freq
    bridge_info_dict['dampingRatio'] = bridge_damp
    bridge_info_dict['modalMass'] = bridge_modal_mass

    file_path = f'{bridge_info_path}.json'
    with open(file_path, 'w', encoding='utf-8') as json_file:
        json.dump(bridge_info_dict, json_file, ensure_ascii=False, indent=4)

    return CAE_df




def vehicle_analysis(param_dict, url_list, CAE_df, KDE_pic_path, \
                     bridge_name, bridge_code):
    Acc_Data = []
    GPS_Data = []

    for url in url_list:
        response = requests.get(url)
        response.raise_for_status()  # 确保请求成功
        content = response.content.decode('utf-8')

        # 假设数据是 CSV 格式，用 pandas 读取
        data = pd.read_csv(io.StringIO(content))

        # 提取数据
        Acc_Data.append(data.iloc[:, 3].to_numpy())
        GPS_Data.append((data.iloc[:, 4].to_numpy(), data.iloc[:, 5].to_numpy()))

    natural_frec = CAE_df['固有频率(Hz)']

    freq_range = [0.8*natural_frec[0],1.2*natural_frec[2]]

    fs = param_dict["采样频率"]
    n_sample = len(url_list)
    band_width = 0.05
    lat_ref1 = GPS_Data[0][0][0]
    lon_ref1 = GPS_Data[0][1][0]
    lat_ref2 = GPS_Data[0][0][-1]
    lon_ref2 = GPS_Data[0][1][-1]

    L_b = param_dict["桥梁全长"]
    N_c = round(L_b / 30)
    EPSG = param_dict["桥梁EPSG代号"]

    fsr = 20
    flo = 0
    fhi = 8
    c0 = 0.3
    alpha = 0.1
    NR = 5
    df = 0.001
    f_inter = np.arange(freq_range[0], freq_range[1] + df, df)
    P_sm_N = np.zeros((len(f_inter), N_c))

    for i in range(n_sample):
        acc_all_z = Acc_Data[i]
        local_x, local_y = gps_to_local(GPS_Data[i][0], GPS_Data[i][1], lat_ref1, lon_ref1, lat_ref2, lon_ref2, EPSG)
        index_ms = (local_x >= 0) & (local_x <= L_b)
        local_x = local_x[index_ms]
        local_y = local_y[index_ms]
        acc_all_z = acc_all_z[index_ms]

        dt = 1 / fsr
        data_proc = dataproc(acc_all_z, fs, fsr, flo, fhi)
        p, q = rat(fsr / fs)
        local_x = local_x[500:-1:q]
        local_y = local_y[500:-1:q]
        t = np.arange(0, len(data_proc)) * dt
        target_length = len(t)
        if len(local_x) > target_length:
            local_x = local_x[:target_length]
        elif len(local_x) < target_length:
            local_x = np.pad(local_x, (local_x[-1], target_length - len(local_x)), mode='constant')
        # 同步压缩小波变换 (SST)
        wavelet = 'morlet'
        # Wx, scales = cwt(data_proc, wavelet)
        Tx, Wx, ssq_freqs, scales = ssq_cwt(data_proc, wavelet=wavelet, fs=fsr)
        Tx = Tx[::-1,:]
        ssq_freqs = ssq_freqs[::-1]
        Tx = Tx[(ssq_freqs >= freq_range[0]) & (ssq_freqs <= freq_range[1]), :]
        ssq_freqs = ssq_freqs[(ssq_freqs >= freq_range[0]) & (ssq_freqs <= freq_range[1])]

        TT, FF = np.meshgrid(t, ssq_freqs)
        TQ, FQ = np.meshgrid(t, f_inter)
        swt_inter = griddata((TT.flatten(), FF.flatten()), np.abs(Tx).flatten(), (TQ, FQ), method='nearest')

        lmax_index = np.zeros_like(swt_inter, dtype=bool)
        for col in range(swt_inter.shape[1]):
            column_data = swt_inter[:, col]
            peaks, _ = find_peaks(column_data)
            lmax_index[peaks, col] = True
    
        logical_matrix = np.zeros_like(swt_inter, dtype=bool)
        for col in range(swt_inter.shape[1]):
            column_data = swt_inter[:, col]
            sorted_indices = np.argsort(column_data)[-round(alpha * len(column_data)): ]
            logical_matrix[sorted_indices, col] = True

        P_lm = np.zeros_like(swt_inter)
        P_lm[lmax_index & logical_matrix] = swt_inter[lmax_index & logical_matrix]      

        P_sm_n = p_sum_r(L_b, c0, N_c, local_x, f_inter, P_lm)
        P_sm_N += P_sm_n

    # plt.imshow(P_sm_N)
    # plt.show()
    # plt.close('all')
    

    # 排序和筛选
    index_fd = np.argsort(P_sm_N, axis=0)[::-1]
    index_fd[NR:, :] = 0
    a = index_fd[index_fd != 0]
    f_s = f_inter[a]
    fs = f_s[(f_s > 0.9 * natural_frec.iloc[0]) & (f_s < 1.1 * natural_frec.iloc[2])]

    # 核密度估计
    kde = gaussian_kde(fs, bw_method=band_width)
    xi = np.linspace(fs.min(), fs.max(), 100)
    f_v = kde(xi)

    threshold = 0
    pks_idx, locs = find_peaks(f_v, prominence=threshold)

    # 匹配频率
    pks_x = xi[pks_idx]

    # ===== 新增绘图代码（20250329-1535） =====
    plt.figure(figsize=(10, 6))
    plt.scatter(pks_x, f_v[pks_idx], c='red', marker='x', label='Peaks')
    plt.plot(xi, f_v, label='KDE', linewidth=2)
    
    plt.title(bridge_code)
    plt.xlabel("Frequency (Hz)")
    plt.ylabel("Density")
    plt.grid(True, alpha=0.3)
    plt.legend()
    
    # 保存图片
    plt.savefig(f'{KDE_pic_path}.png', dpi=300, bbox_inches='tight')
    plt.close('all')  # 关闭图形释放内存
    

    modal_frec = []  # 存储最终结果
    for target in natural_frec:
        distances = np.abs(pks_x - target)
        closest_peak_idx = distances.argmin()  # 找到最近的峰值索引
        modal_frec.append(pks_x[closest_peak_idx])

    for i, modal_res in enumerate(modal_frec):
        if i > 0 and modal_res == modal_frec[i-1]:
            modal_frec[i] = natural_frec[i]
        if abs(modal_res - natural_frec[i]) > 0.03 * natural_frec[i]:
            modal_frec[i] = natural_frec[i]
    
    return modal_frec




if __name__ == "__main__":

    bridge_type = 3
    sample_type = 1
    
    NN_df = pd.DataFrame({
        '振型阶数': [1, 2, 3],
        '固有频率(Hz)': [0.26, 0.33, 0.6],
        '振型阻尼比': [0.01, 0.02, 0.03]
    })

    CAE_df = pd.DataFrame({
        '振型阶数': [1, 2, 3],
        '固有频率(Hz)': [0.26, 0.33, 0.6],
        '振型阻尼比': [0.01, 0.02, 0.03],
        '模态质量(t)': [100, 200, 300]
    })
    

    url_list = [
        "https://haoyaguanjia.oss-cn-hangzhou.aliyuncs.com/bridge/dev/files/csv/1001/1001_sample_02.csv",
        "https://haoyaguanjia.oss-cn-hangzhou.aliyuncs.com/bridge/dev/files/csv/1001/1001_sample_01.csv"
    ]

    param_dict = {
        '采样频率': 100,
        '桥梁全长': 1280,
        '桥梁EPSG代号': 32651,
        'sample_type': sample_type
    }

    res = vehicle_analysis_module(param_dict, url_list, NN_df, CAE_df, bridge_type)
    print(res)
    